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' Abstract. Today, the motion of spacecraft is still described by the classical Newtonian equations 
of motion plus some relativistic corrections. This approach might become cumbersome due to 
the increasing precision required. We use the Relativistic Motion Integrator (RMI) approach to 
numerically integrate the native relativistic equations of motion for a spacecraft. The principle 
of RMI is presented. We compare the results obtained with the RMI method with those from 
the usual Newton plus correction approach for the orbit of the BepiColombo (around Mercury) 
, and Mars-NEXT (around Mars) orbiters. Finally, we present a numerical study of RMI and we 
show that the RMI approach is relevant to study the orbit of spacecraft. 
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1. The Relativistic Motion Integrator (RMI) and an analytical 
5h ! development 

The software RMI presented in Pireaux et al. (2006)| [Pireaux et al. (2008)| numerically 



integrates the relativistic equations of motion 
; dr 2 dr dr 1 ' ' 

in 



for a given metric G^ v where X^ = (cT, X, Y, Z) are the coordinates, r is the proper 
time and Y^ v are the ChristofTel symbols of the metric considered, derived numerically. 
As an example, we use this integrator with the planetocentric metric advised by the IAU 



2000 resolutions, described in Soffel et al. (2003) and characterised by a scalar potential 
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' W and a vector potential W\ Until now, we considered only the central body of mass 
. M, so that we have: 

f W(X a ) = ^ [i + J2Z2EL= (^) 1 Pirn(cos9)(Q m cosmcb + S lm smmq 

X : 1 w\x a ) = -f 

(1.2) 

where G is Newton's gravitational constant, c the speed of light, C; m and Si m are related 
to the central gravity field, R e is the equatorial radius of the central body, while the 
vector S is its spin moment and R = \JX 2 + Y 2 + Z 2 . 

It is possible to develop analytically the equations of motion (jl.ip at first Post- 
Newtonian (1PN) order. Doing so, one gets = —S^R + corr, where the correc- 
tions are composed of different types of forces: a Newtonian correction coming from 
the harmonics (proportional to Ci m or Si m ), a relativistic Schwarzschild acceleration 
(proportional to 1/c 2 ), a relativistic correction coming from the harmonics (propor- 
tional to Cimjc 2 or Si m /c 2 ), a relativistic coupling between harmonics (proportional 
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Figure 1. Corrections due to the relativistic Figure 2. Corrections due to the relativistic 
Schwarzschild acceleration. contribution from the first harmonics (C20, C22 

and #22). 



to Ci m Ci m /c 2 , Ci m Si m /c 2 or Si m Si rn /c 2 ) and finally a relativistic Lense-Thirring accel- 
eration (proportional to the spin momentum over c 2 ). 



2. Results for the BepiColombo and Mars-NEXT missions 

The analytical development described above is used to assess the order of magnitude 
of each separate effect and to validate the RMI method. Figures [TJ [2] and [3] show the 
separate impact of the relativistic effects in terms of cartesian coordinates (X, Y, Z) . 
radial distance R and L = to + w, where lj is the argument of pericenter and w is 
the true anomaly. The orbital parameters of the BepiColombo mission can be found 
in Balog et al. (2000) a = 3389 km, e = 0.162, i = 90 ° . The orbital parameters of the 
Mars-NEXT mission (see |Chicarro et al. (2008)| are: a = 3896 km, e = and i = 75 ° . 
The numerical integration has been performed over 5 orbital periods. 



3. Numerical precision 

The numerical derivative has to be treated carefully. In the implementation of RMI, 
we used a fourth-order numerical derivative 

f(x - 2h) - 8/(x -h) + 8f(x + h)- f(x + 2h) 



D h + 0(h 4 ) 



12h 



+ 0(h 4 ). (3.1) 



For large h, the discretization error is important, while for small h, the roundoff error 
increases. We use an optimal derivation step computed analytically for a function ~ (see 
Figured]), given by h opt = (45ea 6 c 2 /(960GAf ^Z 5 ( |Kincaid and Cheney (2002)] ) where 
e is the machine precision and a is the semi- major axis. As can be seen on Figured we 



derive H 



G 



with 



the Minkowski metric, instead of G^ v since it is more 



stable numerically. Moreover, we use a Richardson extrapolation in order to increase the 
precision on the derivative (Richardson (1927)). This extrapolation uses two estimations 
of the derivative of order 4 with different step size (Dh and with k a real factor) to 
construct an estimation of order 8 



f(x) 



k A Dh/k - Dh 
fc 4 - 1 



0(h 8 ). 



(3.2) 



A Relativistic Motion Integrator (RMI) 



3 



5e-06 
4e-06 
3e-06 
2e-06 
1e-06 
Oe+00 
-1e-06 
-2e-06 
-3e-06 
-4e-06 



BepiColombo 



10000 15000 20000 25000 30000 35000 
Planetocentric Coordinate Time (s) 




derivative of G v 
derivative of H 
Richardson extrapolation 




5000 10000 15000 20000 25000 30000 
Planetocentric Coordinate Time (s) 



Mm) 



Figure 4. Relative precision of the metric 



on 4.- j j. ii w ■ i- derivative as a function of the discretization step 
r igure 3. Corrections due to the relativistic . ^ 

Lense-Thirring acceleration. 



Typically, the value of k is often chosen as 2 or 1.5. We can see on Figure Q] that this 
procedure can increase the derivative precision. It is furthermore less sensitive to the 
choice of /i pt- With such an implementation, it is possible to show that the relative 
precision of RMI is of the order of 10~ 12 in double precision and 10~ 22 in quadruple 
precision. 



4. Conclusion 

We have shown that RMI is useful to compute relativistic orbits for different missions. 
It includes all the relativistic effects (up to the corresponding order of the metric) . It is 
quite easy to use, since the user only has to change the metric module if he wants to 
change the metric. Until now, RMI is only a prototype that is more time consuming than 
the usual 1PN approach. Nevertheless, it is possible to reduce drastically the integra- 
tion time required by the RMI method via proper coding and using parallelization (the 
computation of the Christoffel symbols can easily be parallelized). 
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